Automatic selection of multiple response functions for generalized Richardson-Lucy spherical deconvolution of diffusion MRI data

Tese de mestrado integrado em Engenharia Biomédica e Biofísica (Sinais e Imagens Médicas), Universidade de Lisboa, Faculdade de Ciências, 2021O processo de desenvolvimento do cérebro humano tem sido objeto de estudo desde há vários anos, levando a avanços significativos no que diz respeito à compreensão das suas diferentes fases e mecanismos. Visto que este desenvolvimento resulta de uma série de complexos processos dinâmicos e adaptativos, existe uma busca contínua de informação sobre a sua organização estrutural e funcional, bem como o seu processo de maturação. A ressonância magnética de difusão (dMRI) é uma técnica bastante completa no que diz respeito à análise do cérebro in vivo. Esta técnica é utilizada para realizar um mapeamento quantitativo, através da aplicação de modelos como o modelo de difusão tensorial (DTI). Estes modelos fornecem medidas que caracterizam o cérebro, tais como a anisotropia fraccional (FA) e difusividade média (MD), permitindo assim a quantificação de microestruturas e consequentemente a reconstrução de feixes de substância branca (WM) que ligam diferentes regiões cerebrais. Dadas as suas propriedades de difusão anisotrópica e a sua constituição fibrosa, as fibras de WM têm sido amplamente estudadas através da dMRI. Além disso, a tractografia tornou-se a abordagem padrão no que diz respeito à avaliação da conectividade cerebral usando dados de dMRI. Os métodos de desconvolução esférica (SD) estão entre os mais utilizados para quantificar a distribuição da orientação das fibras (FOD) a partir de dados dMRI do cérebro, sendo que a forma mais comum de o fazer é com desconvolução esférica limitada (CSD). A ideia original da CSD baseia-se no facto de podermos escolher uma função de resposta (RF) representativa de um determinado tecido presente no cérebro e aplicar a SD para resolver o problema de cruzamento de fibras que o modelo de DTI não consegue resolver. Uma vez que o cérebro possui uma complexa organização de tecidos, múltiplos tecidos devem ser considerados. Não é apropriado usar uma RF de WM em todo o cérebro, pois isso pode levar a reconstruções imprecisas da orientação das fibras e a um mau desempenho durante o processo de tractografia. Ao ter em conta múltiplos tecidos, as propriedades da substância cinzenta (GM) e do líquido céfalo-raquidiano (CSF) podem ser quantificadas, e os efeitos de volume parcial (PVE) podem ser reduzidos. Nos últimos anos, tem sido possível adquirir dados “multi-camada” mais complexos e de elevada resolução, mesmo em recém-nascidos, o que permitiu melhorar a técnica de CSD. Consequentemente, esta aquisição também vai melhorar a reconstrução da FOD no cérebro adulto, pois considera os PVE entre diferentes tipos de tecidos. No cérebro neonatal existem algumas diferenças, pois este é constituído por WM em diferentes fases de maturação, e a GM possui características diferentes em comparação com um cérebro adulto. A possibilidade de distinguir diferentes tipos de fibras apenas com base nas suas características microestruturais deve-se às diferenças presentes no cérebro enquanto este se encontra numa fase de desenvolvimento. Em cérebros adultos, é menos provável conseguir observar tais diferenças. Uma das melhores formas de compreender e estudar estes processos de desenvolvimento cerebral é através do estudo do cérebro de neonatais. Como seria de esperar, o cérebro de um recém-nascido não se encontra completamente maturado, sofrendo por isso diversas alterações até estar totalmente desenvolvido. Estas mudanças vão desde o aumento do tamanho do cérebro a alterações ao nível vascular, levando consequentemente a uma alteração dos processos de cognitivos. Em última análise, a aplicação de CSD a dados de “multi-camada” leva a uma extração mais precisa da FOD que por sua vez irá melhorar o processo de tractografia e levará, consequentemente, a uma melhor compreensão do cérebro humano e do seu desenvolvimento, particularmente se aplicada em recém-nascidos e comparada com adultos. O método Generalized Richardson-Lucy (GRL) pode superar os problemas encontrados pela CSD através da realização de SD robusta, suprimindo picos imprecisos na FOD em dados “multi-camada” de dMRI. Este método pode definir múltiplos tecidos que irão aumentar a precisão da estimativa da FOD. No entanto, no método GRL, as três classes de tecidos representadas (WM, GM e CSF) são pré-definidas com valores FA e MD retirados da literatura. Este estudo consistiu em desenvolver um método que determina automaticamente o número de classes (tecidos) necessárias para aplicar corretamente GRL no cérebro com dados “multi-camada”, utilizando para isso os seus valores de FA e MD. O objetivo é aplicar corretamente o método de GRL no cérebro com as classes obtidas, de forma avaliar se existe uma melhoria no processo de estimação das FOD e por sua vez no processo de tractografia. Os dados utilizados neste trabalho consistem em dados de dMRI de dez neonatais e dez adultos, fornecidos pelo Developing Human Connectome Project (dHCP) e pelo Human Connectome Project (HCP), respetivamente. Estes dados já se encontravam num formato pré-processado, pelo que não foi necessário realizar qualquer etapa adicional neste sentido. A primeira parte do estudo consistiu no desenvolvimento do método de deteção automática do número de tipos de tecidos no cérebro. Para isso, todos os dados foram processados no ExploreDTI, um programa de interface gráfica para dados de dMRI e que permite, por exemplo, a realização de tractografia. Este programa foi também usado para extrair os valores de FA e MD dos dados de dMRI dos cérebros dos neonatais e dos adultos, de modo a analisar a sua distribuição de valores por todo o cérebro através de histogramas. De seguida foi aplicado um gaussian mixture model (GMM) aos histogramas de FA e MD, utilizando o MATLAB R2018a, de forma a decompor os dados em classes. Depois de aplicar o GMM aos dados, foi determinado o número ideal de Gaussianas para os mapas de FA e MD. Para isso foi calculado o Bayesian information criterion (BIC) de cada modelo, em que cada um destes se caracteriza por um certo número de Gaussianas. De seguida, foi calculada a probabilidade do valor de cada voxel pertencer a uma das classes escolhidas de FA e MD, atribuiu-se assim uma classe a cada voxel. Posteriormente selecionaram-se as três melhores combinações de FA e MD de cada classe com base na frequência de ocorrência de cada combinação, sendo que cada classe foi definida pela média e desvio padrão das respetivas Gaussianas. Por fim, foram criados mapas espaciais do cérebro com as classes finais, utilizando o MATLAB R2018a. Na segunda parte do estudo aplicou-se o método GRL aos dados, de forma a estimar a RF de cada um dos tecidos que foram selecionados na primeira parte. Estas duas partes do trabalho integram a nossa abordagem, sendo esta designada por "GRL-auto". No método GRL, a RF da GM e do CSF é baseada em valores de FA e MD retirados da literatura, enquanto que o método GRL-auto desenvolvido neste estudo estima esses valores através da seleção automática dos valores de FA e MD que são característicos de cada um destes tecidos. Obtiveram-se os mapas das frações de sinal da WM, GM, e CSF e foram feitas comparações entre o método GRL e GRL-auto. As FOD da WM obtidas com ambos os métodos foram comparadas entre si em regiões de cruzamento de fibras, tanto para neonatais como para os adultos. Por fim, para ambos os métodos, procedeu-se à tractografia em neonatais. Os resultados indicam que, tanto para recém-nascidos como para adultos, existe consistência em relação aos valores de FA e MD e ao seu respetivo número de classes selecionadas. Além disso, conseguem ser observadas diferentes fases de maturação de WM nos neonatais, mas também algumas imperfeições à volta dos ventrículos e regiões onde ocorre cruzamento de fibras. Todos os mapas espaciais de FA e MD fizeram sentido anatomicamente, sendo consistentes quer nos neonatais quer nos adultos, demonstrando assim a eficácia deste método. Os mapas de sinal das frações de WM, GM, e CSF apresentaram valores plausíveis e concordância com a anatomia esperada, para além de consistência tanto nos recém-nascidos como nos adultos. Os mapas de frações de sinal dos adultos praticamente não apresentaram diferenças entre os dois métodos. No entanto, os neonatais mostraram algumas diferenças notáveis, particularmente nos mapas de GM e CSF. Os resultados relativos às FODs não mostraram diferenças significativas no que diz respeito aos adultos. No entanto, para os neonatais, o método GRL-auto estimou FODs de elevada qualidade na WM, em comparação com o método GRL. Além disso, o método GRL-auto detetou mais picos plausíveis em regiões de cruzamento de fibras par além de uma diferença angular maior entre os principais picos das FOD, em comparação com o método GRL. Por fim, este método demonstrou uma melhoria no processo de tractografia, o que por sua vez levará a uma melhor compreensão do cérebro humano e do seu desenvolvimento. Conclui-se assim que o método desenvolvido neste estudo é eficiente e mostra consistência no que diz respeito ao processo de seleção automática do número de tecidos necessários para efetuar CSD no cérebro. Observou-se uma melhoria na tractografia das fibras, o que permitirá uma melhor compreensão da maturação do cérebro bem como das conexões entre as diversas regiões, tendo-se, assim, cumprido o objetivo principal deste trabalho.To understand the development of the human brain, more detailed information is required regarding the structural and functional cerebral organization and maturation. This development is the product of a complex series of dynamic and adaptive processes, and one of the best ways to understand it is through the study of the neonatal brain. The neonatal brain is not fully developed as it would be expected, so it goes through many changes regarding brain size, vasculature, and cognition. Constrained spherical deconvolution (CSD) is a widely used approach to quantify the fiber orientation distribution (FOD) from diffusion magnetic resonance imaging (dMRI) data of the brain, which allows the reconstruction of more complex white matter (WM) bundles in vivo, including in neonates. However, this method estimates the response function (RF) based on the model of a single fiber population and uses it to try to reconstruct the local WM orientations. Since the brain has a complex tissue organization, multiple tissues must be considered. It is not appropriate to use a WM RF throughout the whole brain because this can lead to spurious fiber orientation reconstructions and bad performance during fiber tractography. By accounting for multiple tissues, properties of grey matter (GM) and cerebrospinal fluid (CSF) can be captured, and partial volume effects (PVE) reduced. The acquisition of more comprehensive high-resolution multi-shell dMRI data offers opportunities to take into account multiple tissue types. Ultimately, these improve fiber tractography and consequently lead to a better understanding of the human brain and its development. The generalized Richardson-Lucy (GRL) method can overcome these challenges by performing robust spherical deconvolution (SD) and suppress spurious FOD peaks on multi-shell dMRI data due to PVE. However, in the GRL method, three tissue classes are typically pre-defined to represent WM, GM, and CSF, using fractional anisotropy (FA) and mean diffusivity (MD) values taken from literature. These two metrics are derived from the diffusion tensor model (DTI), with FA measuring how anisotropic is the tensor in each voxel and MD measuring the average of the diffusion rate at each voxel. This study aims to develop a method that automatically determines the number of tissue types (classes) that are needed to properly perform GRL in each analyzed brain dataset. The dataset used in this work consists of ten neonates and ten adults from the Developing Human Connectome Project (dHCP) and the Human Connectome Project (HCP), respectively. The first part of this study consisted of developing a method for the automatic detection of the number of tissue types in the brain, by applying a gaussian mixture model (GMM) and the Bayesian information criterion (BIC) to automatically extract the number of tissue classes from the histogram of dMRI properties. In the second part, the GRL method was applied to the data to estimate the RF of each tissue that was automatically chosen in the first part, and therefore calculate the FOD and perform fiber tractography. This approach was designated by “GRL-auto”. Lastly, a comparison between the basic GRL formulation and GRL-auto was done. Since GRL uses predefined values calibrated on HCP data, it becomes clear that small differences were expected on such dataset, whereas on dHCP larger differences were expected. Our analysis showed that our method automatically identified three classes in the FA histogram and two classes in the MD histogram when using HCP and dHCP data. Therefore, these results demonstrated consistency regarding the FA and MD values and their respective number of selected classes, for both datasets. Furthermore, different stages of WM maturation were detected in the dHCP data, but also some imperfections around the ventricles and crossing fibers areas. All FA and MD spatial maps were in line with anatomical correspondence and were consistent across all neonatal and adult subjects, demonstrating the efficiency of this method. The values of the WM, GM, and CSF fraction maps were plausible, in line with the expected anatomy, and looked consistent on both HCP and dHCP datasets. The signal fraction maps determined with the HCP data showed almost no difference between GRL and GRL-auto. However, in the dHCP data, we observed notable differences, particularly in the GM and CSF maps. Regarding the FOD estimation, our results showed no difference in the HCP data. Nevertheless, for the dHCP data, GRL-auto estimated high-quality FODs in WM, and detected more peaks in crossing fiber regions and a bigger angular difference between the main FOD peaks, as compared to GRL. Lastly, we showed that GRL-auto led to improvements in fiber tractography, which will likely support gaining a better understanding of the human brain and its development. Therefore, we can conclude that the method developed in this study is efficient and consistent in the automatic selection of the number of tissues needed to properly perform GRL in a brain, given multi-shell data, which was the main goal